% Analysis Reduction Error 
clc;
clear all;
close all;
%% function analysis_red_err_alpha(lower,upper,m,n,r,k, alpha)
m = 5; n = 5; r = 5; k = 2;
lower = 0;
upper = 0.9999999999;

step_a = 2^(-m);
step_b = 2^(-n);
step_c = 2^(-r-k);
max_const = ccm_rd(m, n, r, k, 1);
for t = 1:k
    max_const = max_const + 2^(- r - t);
end
%num = (round((upper-lower)/step_a) + 1)*(round((upper-lower)/step_b) + 1);

alpha = 0:0.01:1;
L = length(alpha);
const = 0:step_c:max_const;
K = length(const);
min_max_err = zeros(1, K);

%% Find min_max_error
disp('Find min max');
for i = 1:K %for each const.
    %find the worst case
    fprintf('Const %d of %d\n', i, K);
    i_const = const(i);    
    index = 0;
    max_err = zeros(1, L);
    for j = 1:L%for each prior
        fprintf('Prior %d of %d \n', j, L);
        j_alpha = alpha(j);        
        for a = lower:step_a:upper
            for b = lower:step_b:upper
                %index = index + 1
                p_ab = prob_joint(a, b, m, n, j_alpha);%probability of a,b        
                %answers true and with const.
                ans_true = a * b;
                ans_const = trunc_mult(a, b, i_const, m, n, r, k);
                %errors
                err_const = ans_const - ans_true;
                %average errors
                max_err(j) = max_err(j) + err_const * p_ab;                
            end
        end
        max_err = abs(max_err);
    end
    min_max_err(i) = max(max_err);
end


%%
figure(1);
ulp = 2^(-r);
plot(const, min_max_err./ulp); 
title('Minimax Total Error (m = 5, n = 5, r = 5, k = 2)');
xlabel('Correction Constant');
ylabel('Total Error(ulp)');
%% Find max_min_error
max_min_err = zeros(1, L);
for i = 1:L %for each prior.
    %find the best solution
    fprintf('Prior %d of %d \n', i, L);
    i_alpha = alpha(i);    
    min_err = zeros(1, K);
    for j = 1:K%for each const
        fprintf('Const %d of %d \n', j, K);
        j_const = const(j);        
        for a = lower:step_a:upper
            for b = lower:step_b:upper
                %index = index + 1
                p_ab = prob_joint(a, b, m, n, i_alpha);%probability of a,b        
                %answers true and with const.
                ans_true = a * b;
                ans_const = trunc_mult(a, b, j_const, m, n, r, k);
                %errors
                err_const = ans_const - ans_true;
                %average errors
                min_err(j) = min_err(j) + err_const * p_ab;                
            end
        end
        min_err = abs(min_err);
    end
    max_min_err(i) = min(min_err);
end
%% Display results
figure(1);
plot(const, min_max_err); 
figure(2);
plot(alpha, max_min_err);

